

clear
/*#delimit;*/
set mem 1000m
pause on
set more off

cd "FILEPATH"

local pathfile = "FILEPATH"

********************************************************************************

use "`pathfile'/DataForJMP_Matching.dta"

gen missing_bonus = 1 if bonusfull == 0 | bonusfull == .
replace missing_bonus = 0 if missing_bonus == .

ttest fulllease, by(missing_bonus) unequal
ttest royalty, by(missing_bonus) unequal
ttest externlease, by(missing_bonus) unequal
ttest environ_bun, by(missing_bonus) unequal
ttest termmonths, by(missing_bonus) unequal

tab instyear missing_bonus

reg royalty missing_bonus

bysort firm_1 wellpadID: egen firmsignedWP = sum(counter) if wellpadID != .
bysort wellpadID: egen totalsignedWP = sum(counter) if wellpadID != .
gen firmshareWP  = firmsignedWP/totalsignedWP

bysort firm_1 api: egen firmsignedApi = sum(counter) if api != ""
bysort api: egen totalsignedApi = sum(counter) if api != ""
gen firmshareApi = firmsignedApi/totalsignedApi

bysort wellpadID firm_1: gen firmWP_ct = _n if wellpadID != .
bysort api firm_1: gen firmApi_ct = _n if api != ""
bysort api: gen api_ct = _n if api != ""

********************************************************************************
* Table 1: Lease Clause and Bundle Summary Statistics

local legal   = "forcemaj pugh verticalpugh insurindem_bun totreport_bun postprod_bun lesseeDefends"
local bads    = "subsurfease_bun nolitigationlessor freewater"
local surface = "surfdamage_bun nosurf_bun surfrest_bun setback"
local extern  = "environ_bun noise_bun freshwaterprotect surfacecasing compressionstation"
local water   = "freshwaterprotect waterquality watersediment nowateruse groundwater_bun"

sum royalty termmonth bonusfull bonus_acre daysaftercess
corr fulllease royalty termmonth bonus_acre bonusfull
corr fulllease royalty termmonth bonusfull badslease externlease surfacelease legallease waterlease

sum fulllease badslease `bads' externlease `extern' surfacelease `surface' legallease `legal' waterlease `water'

********************************************************************************
* Table 2: Parcel Summary Statistics

sum land_sqft abst_no_app_med near_dist_pipe well_lateral_dist parcelIntLat

********************************************************************************
* Table 3: Firm and Pricing Summary Statistics

sum near_dist_pipe well_lateral_dist

bysort firm_1: gen firm_ct = _n
replace operatorLarge = 0 if operatorLarge == .
sum landman operatorLarge if firm_ct == 1
sum landman operatorLarge

bysort instyear instmonth instday: gen year_ct = _n
sum px_last_wti18 px_last_ng18 if year_ct == 1

********************************************************************************
* Table 4: Summary Statistics: Instrumental Variables Mode

replace parcelToAPI = parcelToAPI/1000
sum firm_clust_share_time_cumCMon if royalty != .
sum op_proddir_wells_2003 op_prodgas_wells_2003 op_prod_wells_2003 prodByAbst03_near if royalty != .
sum land_sqft abst_no_app_med parcelToAPI if royalty != .

********************************************************************************
* Table 5, 6, 7: Lease Quality Estimates - Auxiliary Terms

xtset wellpadID

local inst     = "firm_clust_share_time_cumCMon"
local infrast  = "op_proddir_wells_2003 prodByAbst03_near"

local price   = "px_last_wti18 px_last_ng18 pred_wti_vol18 pred_ng_vol18"
local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name `price'"

foreach leasevar of varlist royalty royalty_land termmonths postprodcost {
xi: xtreg `leasevar' `inst' `othercontrol',  fe 
xi: xtreg `inst' `infrast' `othercontrol' if `leasevar' != .,  fe 
xi: xtreg `leasevar' `infrast' `othercontrol', fe 
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first cluster(wellpadID) gmm2s 
}

local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name"

foreach leasevar of varlist fulllease legallease surfacelease badslease externlease daysaftercess {
xi: xtreg `leasevar' `inst' `othercontrol',  fe 
xi: xtreg `inst' `infrast' `othercontrol' if `leasevar' != .,  fe 
xi: xtreg `leasevar' `infrast' `othercontrol', fe 
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first cluster(wellpadID) gmm2s 
}

local inst = "firmcount_wp_Mcum"

local leasevar = "royalty"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s
local leasevar = "bonus_acre"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s
local leasevar = "fulllease"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s
local leasevar = "legallease"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s
local leasevar = "badslease"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s
local leasevar = "externlease"
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', first gmm2s


********************************************************************************
* Appendix Tables 5, 6, 7: Lease Quality Estimates - Auxiliary Terms

local inst     = "firm_abst_share_time_cumAMon"
local infrast  = "op_proddir_wells_2003 prodByAbst03_near"
local price   = "px_last_wti18 px_last_ng18 pred_wti_vol18 pred_ng_vol18"
local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name `price'"

foreach leasevar of varlist royalty royalty_land termmonths postprodcost {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
}
local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name"
foreach leasevar of varlist fulllease legallease surfacelease badslease externlease daysaftercess {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
}

local inst     = "hhi_clust_time_cumMon"
local infrast  = "op_proddir_wells_2003 prodByAbst03_near"
local price   = "px_last_wti18 px_last_ng18 pred_wti_vol18 pred_ng_vol18"
local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name `price'"

foreach leasevar of varlist royalty royalty_land termmonths postprodcost {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
}

local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name"
foreach leasevar of varlist fulllease legallease surfacelease badslease externlease daysaftercess {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol', first  fe 
}

********************************************************************************
* Figure 2: Bundle versus Clause Estimates - Probit IV GMM

xtset wellpadID

local xt = 1
local iv = 1
local IV = "-IV"
local FE = "-FE"
local modeltype = "xtiv_probit"
local bundletype = "full"

local full  = "environ_bun noise_bun freshwaterprotect surfacecasing compressionstation traffic_bun surfdamage_bun nosurf_bun surfrest_bun setback forcemaj pugh verticalpugh insurindem_bun totreport_bun postprod_bun lesseeDefends subsurfease_bun nolitigationlessor freewater"
local fulllabel = `" 1 "Envronmental" 2 "Noise" 3 "Freshwater" 4 "Surface_Casing" 5 "Compression_Station" 6 "Traffic" 7 "Surface_Damage" 8 "No_Surface_Access" 9 "Surface_Restriction" 10 "Setback" 11 "Force_Majuere" 12 "Pugh" 13 "Vertical_Pugh" 14 "Insurance" 15 "Total_Report" 16 "Post_Production_Cost" 17 "Firm_Defends_Title" 18 "Subsurface_Ease" 19 "No_Litigation" 20 "Free_Water" "'

local extern  = "environ_bun noise_bun freshwaterprotect surfacecasing compressionstation traffic_bun"
local externlabel = `" 1 "Envronmental" 2 "Noise" 3 "Freshwater" 4 "Surface_Casing" 5 "Compression_Station" 6 "Traffic" "'

local surface = "surfdamage_bun nosurf_bun surfrest_bun setback"
local surfacelabel = `" 1 "Surface_Damage" 2 "No_Surface_Access" 3 "Surface_Restriction" 4 "Setback" "'

local legal   = "forcemaj pugh verticalpugh insurindem_bun totreport_bun postprod_bun lesseeDefends"
local legallabel  = `" 1 "Force_Majuere" 2 "Pugh" 3 "Vertical_Pugh" 4 "Insurance" 5 "Total_Report" 6 "Post_Production_Cost" 7 "Firm_Defends_Title" "'

local water   = "freshwaterprotect watersediment nowateruse groundwater_bun "
local waterlabel   = `" 1 "Freshwater" 2 "Water_Sediment" 3 "No_Water" 4 "Groundwater" "'

local bads    = "subsurfease_bun nolitigationlessor freewater"
local badslabel   = `" 1 "Subsurface_Ease" 2 "No_Litigation" 3 "Free_Water""'

local othercontrol = "land_sqft parcelIntLat abst_no_app_med parcelToAPI i.instyear i.city_name"

local inst     = "firm_clust_share_time_cumCMon"
local infrast  = "op_proddir_wells_2003 prodByAbst03_near"

local leasevar = `"`bundletype'lease"'

gen plotcount = _n

if `xt' == 1 & `iv' == 0 {
xi: xtreg `leasevar' `inst' `othercontrol' if `leasevar' != . , fe
}
if `xt' == 0 & `iv' == 1 {
xi: ivreg2 `leasevar' (`inst' = `infrast') `othercontrol', vce(cluster wellpadID) gmm2s
}
if `xt' == 1 & `iv' == 1 {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol' , fe
}
else {
xi: reg `leasevar' `inst' `othercontrol' if `leasevar' != ., vce(cluster wellpadID) 
}

qui: matrix list e(b)
qui: matrix list e(V)

gen coef = _b[firm_clust_share_time_cumCMon]
gen stdev = _se[firm_clust_share_time_cumCMon]
gen lower = coef - 1.96*stdev
gen upper = coef + 1.96*stdev

gen coef_clause = .
gen stdev_clause = .
local i = 1

foreach leasevar of varlist ``bundletype'' {
if `xt' == 1 & `iv' == 0 {
xi: xtreg `leasevar' `inst' `othercontrol' if `leasevar' != . , fe
}
if `xt' == 0 & `iv' == 1 {
xi: ivprobit `leasevar' (`inst' = `infrast') `othercontrol', vce(cluster wellpadID)
}
if `xt' == 1 & `iv' == 1 {
xi: xtivreg `leasevar' (`inst' = `infrast') `othercontrol' , fe
}
else {
xi: probit `leasevar' `inst' `othercontrol' if `leasevar' != ., vce(cluster wellpadID) 
}

qui: matrix list e(b)
qui: matrix list e(V)
replace coef_clause = _b[firm_clust_share_time_cumCMon] if plotcount == `i'
replace stdev_clause = _se[firm_clust_share_time_cumCMon] if plotcount == `i'
local i = `i' + 1
}
gen lower_clause = coef_clause - 1.96*stdev_clause
gen upper_clause = coef_clause + 1.96*stdev_clause

gen yline = 0

forval j = 1(1)`i' {
local obs = `i' + `j'
replace plotcount = `j' + 0.5 in `obs'
}
sort plotcount

graph twoway (line yline plotcount if plotcount < `i' ,  lcolor(red)) /*
*/  (rarea lower_clause upper_clause plotcount if plotcount < `i' ,  lcolor(black) lwidth(vthin) color(gray) fin(30)) /*
*/  (connected coef_clause plotcount if plotcount < `i' ,  lcolor(black) mcolor(black) msize(vsmall) ) /*
*/  (rarea lower upper plotcount if plotcount < `i' ,  lcolor(black) lwidth(vthin) color(none) fin(30)) /*
*/  (connected coef plotcount if plotcount < `i' ,  lcolor(black) mcolor(black) msize(vsmall) ) /*
*/ , graphregion(color(white)) bgcolor(white) ylabel(,nogrid) legend(order(2 "Clause" 4 "Bundle") size(vsmall)) /*
*/	ytitle("Probit`IV'`FE' estimates, Concentration", size(small)) xtitle("", size(small)) xlabel(``bundletype'label', angle(45))   /*
*/  name(hold, replace) yline(0) /*
*/  ylabel(,labsize(vsmall)) xlabel(, labsize(vsmall))

drop lower* upper* plotcount coef* stdev* yline

graph combine hold, saving(hold,replace)
graph export `bundletype'`modeltype'.png, replace

********************************************************************************
* Table 8: Geographic Market Characteristics

gsort api -firmshareApi
by api: gen highshare = _n
sum firmshareApi if highshare == 1
sum totalsignedApi if highshare == 1

drop api *Api highshare
duplicates drop

gsort wellpadID -firmshareWP
by wellpadID: gen highshare = _n
sum firmshareWP if highshare == 1
sum totalsignedWP if highshare == 1

bysort wellpadID: egen land_sqft_min = min(land_sqft)
bysort wellpadID: egen land_sqft_max = max(land_sqft)

gen land_sqft_pct = land_sqft_max/land_sqft_min - 1
sum land_sqft_pct if wellpad_ct == 1, det

sum totalsignedWP, det
gen well_less_50  = 1 if totalsignedWP < r(p50)
replace well_less_50 = 2 if well_less_50 == .
gen well_less_10  = 1 if totalsignedWP < r(p10)
replace well_less_10 = 2 if totalsignedWP > r(p90)
gen well_less_25  = 1 if totalsignedWP < r(p25)
replace well_less_25 = 2 if totalsignedWP > r(p75)

* Change to 50, 25, 10 generate the Market Heterogeneity results for each column
local pct = 50
local first = 1
foreach var of varlist land_sqft well_dist lateral_dist parcelIntLat abst_no_app_med fulllease {

sum `var' if well_less_`pct' == 1
scalar lowerval = r(mean)
sum `var' if well_less_`pct' == 2
scalar upperval = r(mean)
scalar diff = (upperval/lowerval-1)*100
sum `var' if well_less_`pct' == 2
matrix tempmat = r(N)\lowerval\upperval\diff

if `first' == 1 {
matrix tempfull = tempmat'
local first = 2
}
else {
matrix tempfull = tempfull\tempmat'
}
}
matrix list tempfull

sum firmsigned, det
gen firm_less_50  = 1 if firmsigned < r(p50)
replace firm_less_50 = 2 if firm_less_50 == .
gen firm_less_10  = 1 if firmsigned < r(p10)
replace firm_less_10 = 2 if firmsigned > r(p90)
gen firm_less_25  = 1 if firmsigned < r(p25)
replace firm_less_25 = 2 if firmsigned > r(p75)

* Change to 50, 25, 10 generate the Firm Heterogeneity results for each column
local pct = 10
local first = 1
foreach var of varlist land_sqft well_dist lateral_dist parcelIntLat abst_no_app_med fulllease {

sum `var' if firm_less_`pct' == 1
scalar lowerval = r(mean)
sum `var' if firm_less_`pct' == 2
scalar upperval = r(mean)
scalar diff = (upperval/lowerval-1)*100
sum `var' if firm_less_`pct' == 2
matrix tempmat = r(N)\lowerval\upperval\diff

if `first' == 1 {
matrix tempfull = tempmat'
local first = 2
}
else {
matrix tempfull = tempfull\tempmat'
}
}
matrix list tempfull

********************************************************************************
* Appendix Table 11: Bonus "Production Function"

bysort abst_no: egen land_avg = mean(land_sqft)
 
qui: tab city_name, gen(city_num)
qui: tab firm_1, gen(firm_num)
qui: tab fdate_inst, gen(fdate_inst)

gen share_land_size = firmshare_clust*contig_acres 

local controls = "land_sqft parcelIntLat abst_no_app_med land_avg"

xtset wellpadID
xi: xtreg bonusfull `controls' i.instyear firm_num1-firm_num20 firm_num22-firm_num27, fe
predict bonus_hat, xb

sum bonus_hat, det
replace bonus_hat = . if bonus_hat < 200
replace bonus_hat = . if bonus_hat > 25000
replace bonus_hat = bonus_hat*land_acres

gen exp_royalty_bonus = 0.25*exp_royalty*0.9 + bonus_hat if bonus_hat != .
gen bonus_frac = bonus_hat/exp_royalty_bonus if bonus_hat > 0
